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Abstract 

We discuss a sampling algorithm which generates flat histogram in 
energy. In combination with transition matrix Monte Carlo, the density 
of states and derived quantities such as entropy and free energy as a 
function of temperature can be computed in a single simulation. 



1 Introduction 

The basic problem in equilibrium statistical mechanics is to compute the canon- 
ical average 

_ E M ^)ex P (-g(q)/fcT) 
{)T E W exp(-i%)/fcT) ■ W 
In addition, free energy 

F = -kThi^2exp(-H(a)/kT) (2) 
M 

and related quantity entropy are also very important. Standard Monte Carlo 
method, e.g., Metropolis importance sampling algorithm, is simple and general. 
However, computation of free energy Q is difficult with such method. 

Over the last decade, there have been a number of methods addressing this 
problem g, [|, H [| [|. A common theme in these approaches is to evaluate 
the density of states n(E) directly. If this can be done with sufficient accuracy, 
then the summation over all configuration states can be rewritten as sum over 
energy only, e.g., 

F = -kTln^2n(E)exp(-E/kT). (3) 

E 

Can we evaluate n(E) with a uniform relative accuracy for all El Our experience 
with flat histogram method suggests a positive yes. 
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2 Broad histogram equation 



Oliveira et al showed the validity of following equation relating density of 
states with the microcanonical average number of moves in a single spin flip 
dynamics: 

n(E){N(a, AE)) E = n(E + AE)(N(a', -AE)) E+AE . (4) 

This equation is equivalent to a detail balance condition in an infinite temper- 
ature transition matrix Monte Carlo simulation 0. The quantity N(a,AE) is 
the number of ways that the system goes to a state with energy E + AE by 
a single spin flip, given that the current state is er with energy E. The aver- 
age (■ ■ -}e is performed over all the states with a fixed initial energy E (i.e., a 
microcanonical average). 

3 A flat histogram dynamics 

Consider the follow Monte Carlo dynamics 

1. Pick a site at random. 

2. Flip the spin with probability r(E'\E). 

3. Sample N(a,AE), i.e., accumulate the statistics for (N(a, AE)). 

4. Go to 1. 

The flip probability r is given by 

V 1 ' V (N(a,AE)) E J K ' 

where the current state a has energy E, and new state a 1 with one spin flipped 
has energy E' = E + AE. 

With the above choice of flip rate, we can show that detail balance is satisfied 

r{E'\E)P{a) = r(E\E')P(a') (6) 

if -P(cr) = const jn(E (<t)), since this equation is equivalent to the broad his- 
togram equation, Eq. (||). The histogram in energy is then H(E) cx ^2 -P( a ) = 
constant. It turns out that the choice of flip rate is not unique; many other 
formulas are possible 

Due to Eq. ([|), the flip rate is also equal to min[l, n(E)/n(E')]. This is 
exactly Lee's method [Q of entropy sampling (which is equivalent to multi- 
canonical method However, since neither n(E) nor (N(a, AE)) E is known 
before simulation, the way by which simulation gets boot-trapped is quite differ- 
ent. Our method is very efficient in this respect. Another important difference 
is that we take (N(a, AE)) E as our primary statistics, from which we derive 
the density of states n(E). Apart from the number of iterations needed, the 
transition matrix results are in general more accurate than that obtained using 
the energy histogram M. 
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4 Simulation procedures 

The first approximate scheme we use is to replace the true microcanonical av- 
erage by a cumulative sample average 



where the samples a 1 are configurations generated during simulation with energy 
E. Each state in the simulation contributes to N{- ■ ■) for some E. For those 
E' where data are not available, we set r = 1. This biases towards unvisited 
energy states. This dynamics does not satisfy detail balance exactly since the 
transition rate is fluctuating. However, test on small systems shows that the 
sample averages do converge to the exact value with errors proportional to the 
inverse of square root of Monte Carlo steps. 

A two-stage simulation will guarantee detail balance. Stage one is the same 
as described above. In stage two, we adjust the approximate transition matrix 
obtained in stage one such that detail balance is satisfied exactly. In this stage, 
the flip rate is fixed and will not fluctuate with the simulation. The second 
stage simulation is dynamically equivalent to Berg's multicanonical or Lee's 
entropy sampling dynamics. The stage two can be iterated so that the simulated 
ensemble approaches the ideal multicanonical ensemble, but we found two-stage 
or even a single stage simulation already gives excellent results. 

The simulation can also be combined with the N-fold way (T(J with little 
overhead in computer time, since the quantity N(a, AE) needed in the N-fold 
way is already computed. In addition, not only we can have equal histogram 
(multicanonical) , we can also generate "equal-hit" ensemble, where each energy 
of distinct states is visited equally likely 

5 Density of states from transition matrix 

The density of states is related to the transition matrix Te : ae = (N(a, AE)) e /N 
by Eq. where N is the total number of possible moves. Since there are more 
equations than unknown n(E), we use a least-square method to obtain "opti- 
mal" solution. Let S(E) — lnn(E), we consider 



subject to all the conditions known. For example, for the Ising model, we have 
n{E m i n ) = n(E max ) — 2, and n(E) = 2 . The variance a 1 is the variance 
of the quantity ItiTe^ae/Te 1 ,-ae obtained from sets of Monte Carlo data. It 
is also possible to work with the matrix T directly with the conditions that T 
must satisfy. 
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6 Conclusion 



We proposed an algorithm which samples energy E uniformly. Comparing to 
multicanonical simulation, the method offers a very easy way of starting the 
simulation. The dynamic characteristics are similar to well-converged multi- 
canonical Monte Carlo dynamics. For example, the tunneling time for 10-state 
Potts model in two dimensions is about t oc L 2 6 for L x L system. It is very 
easy to combine statistics from several simulations, including parallel simula- 
tions. It is an efficient method for computing density of states and all the related 
thermodynamic quantities. 
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